# -*- coding: utf-8 -*-
from __future__ import unicode_literals

"""
Spyder Editor

This is a temporary script file.
非制造业 生存曲线分析
"""
import pylab as pl
import numpy as np
from scipy.stats import gamma


def show(l_date, k, r):
    x = pl.frange(1, l_date, 1)  # x轴
    z = gamma.pdf(x, k, 0, r)
    y = ((k - 1) / x - 1 / r) * gamma.pdf(x, k, 0, r)  # y轴
    l1 = 0  # 分界点索引
    l2 = np.argmin(y)
    l3 = 0
    r = abs(np.min(y))
    for i in range(l_date - 1):
        if (y[i] * y[i + 1] < 0):
            i1 = abs(y[i + 1] - y[i])
            i2 = abs(y[i + 1])
            l = i2 / i1
            l1 = x[i] * l + x[i + 1] * (1 - l)
            break
    distence = abs(l2 - l1) / 2.0
    for i in range(l_date - 1):
        if i > l2 + distence:
            if (abs(y[i]) + abs(y[i] - y[i + 1])) <= r:
                l3 = i
                break
    if l3 < l2:
        l3 = l_date + 1
    else:
        l3 = x[l3]
    x = map(lambda x:x-1, x)
    return [x, z], [l1, x[l2], l3]

# if __name__ == '__main__':
#     k = 1.42614631205
#     r = 3.95550588726
#     print show(k, r)
